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ABSTRACT 

We analyze the statistical properties of bubble models for the large-scale distribu- 
tion of galaxies. To this aim, we realize static simulations, in which galaxies are 
mostly randomly arranged in the regions surrounding bubbles. As a first test, we 
realize simulations of the Lick map, by suitably projecting the three-dimensional 
simulations. In this way, we are able to safely compare the angular correlation 
function implied by a bubbly geometry to that of the APM sample. Quite remark- 
ably, we find that several bubble models provide an adequate amount of large-scale 
correlation, which nicely fits that of APM galaxies. Further, we apply the statistics 
of the count-in-cell moments to the three-dimensional distribution and compare 
them with available observational data on variance, skewness and kurtosis. Based 
on our purely geometrical constructions, we find that a well defined hierarchical 
scaling of higher order moments up to scales ~ 70/i~^Mpc. We show that this must 
be expected for any non-Gaussian distribution in the weak-coherence regime. The 
overall emerging picture is that the bubbly geometry is well suited to reproduce 
several aspects of large-scale clustering. Furthermore, the statistical test we apply 
are able to discriminate between different models. We find that models with fixed 
bubble radius have troubles to account for the angular correlation of APM galaxies, 
while a shallow spectrum fails to reproduce the observed three-dimensional skew- 
ness. A model with a rather steep spectrum of bubble of radii is rather adequate 
to account for the observational tests that we considered. 
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1 Introduction 

The most striking feature of the large-scale structure of the Universe is the remarkable variety 
of coherent structures and the extension of the involved scales. Its interpretation really rep- 
resents a difficult task for any theoretical model about the origin of primordial fluctuations. 
The most eminent victim of such a largely structured Universe is probably the standard cold 
dark matter (CDM) scenario. In its minimal version, this model assumes initial adiabatic and 
Gaussian fluctuations, dominated by highly non-relativistic non-baryonic particles in a ^2 = 1 
Universe, with galaxies substantially more clustered than the underlying dark matter (see, 
e.g., Blumenthal et al. 1984; see also Davis et al. 1992 and Liddle & Lyth 1993, for recent 
reviews about CDM). Although successful at small and intermediate scales, this model has 
been shown to have serious troubles in accounting for clustering data at scales > 20/i~^Mpc, 
such as clustering of rich clusters ( e.g.. White et al. 1987), large-scale galaxy angular correla- 
tion (Maddox et al. 1990; Collins et al. 1991), large-scale galaxy peculiar motions (Vittorio, 
Juszkiewicz Sz Davis 1986; Bertschinger et al. 1990), count-in-cell statistics of IRAS-QDOT 
(Saunders et al. 1991) and Stromlo-APM (Loveday et al. 1992) galaxies, large-scale power 
spectra (Peacock 1991; Feldman, Kaiser & Peacock 1993), abundance of rich clusters (White, 
Efstathiou & Frenk 1993). 

Several adjustments of the standard CDM model have been proposed in order to overcome 
its deflcit of large-scale power, which can be divided into two main categories; those which 
modify the scheme for the production of primordial inhomogeneities and those which try to 
suitably modify the prescription for galaxy identiflcation out of dark matter fluctuations. The 
flrst category includes models with primordial power-spectrum tilted away from the Zel'dovich 
one, predicted by standard inflation ( e.g., Lucchin & Matarrese 1985; Cen et al.l992), low- 
density Universe with non- vanishing cosmological constant iyimatter — 0.2, ^matter + — Ij 
e.g., Efstathiou, Sutherland & Maddox 1990), non-Gaussian initial fluctuations (Moscardini 
et al. 1991; Weinberg & Cole 1992), inclusion of a hot particle component supplying large- 
scale power (Valdarnini & Bonometto 1985; Achilli, Occhionero & Scaramella 1985; Klypin 
et al. 1992). The second class of non-standard CDM models includes low-bias mechanisms 
of galaxy formation, in which the present time corresponds to a dynamically more evolved 
stage (Couchman & Carlberg 1992; see, however, Moscardini et al. 1993), and "cooperative" 
non-local biassing (Bower et al. 1993), in which galaxy formation is triggered by astrophysical 
mechanisms, able to introduce large-scale coherence ( e.g., inhibition of galaxy formation due 
to high-redshift QSO reionization; Babul & White 1991). 

In this paper we follow a purely geometrical approach by modelling the large-scale clus- 
tering by means of a distribution of spherical bubbles, with galaxies mostly arranged in the 
volume not occupied by such bubbles, while leaving them almost devoid. Therefore, we are not 
concerned with the dynamics underlying the above models for large-scale structure formation. 
Instead, we attempt to parametrize their dynamical content through a suitable geometrical 
construction. 

Bubble models have been already considered in several contexts. Generation of bubbles 
can be provided by invoking primordial phase-transitions, occurred along with the inflationary 
expansion ( e.g., La 1991). The only theory to date for primeval bubble generation, which also 
predicts a deflnite spectrum of bubble sizes, is the extended inflation (La & Steinhardt 1989; 
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for a review see Kolb 1991). In this scenario, bubbles are generated as true- vacuum regions, 
tunnelling out from the initial false- vacuum state and subsequently stretched to cosmological 
sizes by the accelerated expansion. In the original model, the number of bubbles Nb{R) with 
radius larger than R is roughly a power-law. 



where p is related to the microphysical parameters of the phase transition and Rm is the 
maximum bubble radius inside the horizon, also fixed by the model. If the inflation is driven 
by the Brans-Dicke field coupled through the constant 6, like in La & Steinhardt (1989), then 
p — 3 + 4:/{b+l/2). The determination of the slope p of the present bubble spectrum may then 
be related to the fundamental constant b. A potential problem in this mechanism of bubble 
production lies in the predicted Rm values, which, in all the models considered so far, are 
either in conflict with the microwave background isotropy or too small to be of cosmological 
interest (Liddle & Wands 1991). A primordial generation of bubbles clearly fixes non-Gaussian 
initial conditions, since the expectation of finding underdense regions is larger than that of 
finding overdense regions, and the resulting probability density function of the density field is 
negatively skewed. 

A bubbly galaxy distribution naturally arises also in the explosion scenario (Ostriker & 
Cowie 1981; Ikeuchi 1981). Differently from the gravitational instability picture, in this sce- 
nario energy perturbations of non-gravitational origin drive material away from the seeds of 
the explosions, sweeping primordial gas into dense, expanding shells. As these shells cool, 
their fragmentation could give a further generation of objects, which again explode, thus am- 
plifying the process and giving rise to large-scale structure formation. At the end, galaxies 
arc arranged on spherical shells and rich clusters are expected to be placed at the intersection 
of three expanding bubbles. A variety of physical mechanisms may generate such explosions, 
such as supermassive stars or supernovae from the earliest galaxies. However, it is at present 
not clear whether such energy sources are sufficient to create the large voids that are observed. 
Moreover, suitable initial conditions are anyway needed to generate primordial objects, which 
act as seeds of explosions. It is also to be remarked that the recent analysis of the COBE FI- 
RAS spectrum (Wright et al.l993) seems to rule out the explosions scenario, since the energy 
release would have distorted the cosmic blackbody spectrum above the observational bounds. 

A more conventional mechanism for nearly spherical bubble generation is provided by 
gravitational instability in the mildly non-linear regime. In the pancake scenario, overdense 
regions tend to increase their asphcricity during the gravitational collapse, while underdense 
regions expand and become more and more spherical. As a result, the most part of space turns 
out to be occupied by almost devoid spherical regions, with galaxies arranged in sheet-like 
structures around them. It is however clear that in order for gravitational dynamics to produce 
as large voids as observed, a primordial fluctuation spectrum is required, which produces an 
adequate amount of large-scale fluctuation power, probably more than supplied by standard 



Purely geometrical random bubble distributions have been already shown to be able to 
account for some features of large-scale clustering. In the framework of the explosion model, 
Weinberg, Ostriker & Dekel (1989) used a random distribution of expanding shells with both 
a power-law and a uniform distribution of radii. After identifying clusters at the intersection 
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of three of such shells, the statistics of their observed distribution is quite well reproduced. 
A detailed description of the galaxy distribution in terms of Voronoi tessellation (van de 
Weygaert & Icke 1989; Yoshioka & Ikeuchi 1989; van de Weygaert 1993) also reveals that a 
topology dominated by big voids surrounded by a sheet-like galaxy distribution gives a fair 
representation of several clustering features. 

In this paper, we analyze statistical properties for simulations of bubbly galaxy distri- 
butions and compare them with similar tests, already presented in the literature, about the 
observed galaxy clustering. This allows us to check in detail the ability of this purely geo- 
metrical model to account for the large-scale texture of the Universe. In the simulations, the 
centers of spherical bubbles are randomly distributed and the radii are chosen either fixed or 
taken from a power law spectrum, like that of Eq.(|l]), with a lower cut-off radius Rm of the 
order of a few megaparsecs. A list of the parameters involved in each of the models that we will 
consider is presented in Table I. A letter F marks the models with fixed radius size. The model 
parameters are always tuned in such a way that the resulting bubble volume is 90% of the 
total volume. Galaxies are randomly arranged in the space not occupied by bubbles, while an 
unclustered component, a fraction /„ of the total number of particles, is also added to account 
for a residual galaxy population inside the underdense regions. The simulation box bound- 
ary is taken to be periodic. By varying the model parameters, we generate simulated galaxy 
distributions having different amounts of large-scale clustering. It is clear that using such 
simple simulations we completely neglect the presence of gravitational clustering at the scales 
of non-linearity. For this reason, we restrict our analysis to larger scales (> S/i^^Mpc), where 
non-linear gravitational dynamics does not work and the clustering should be accounted for 
just by the large-scale geometry. The remarkable advantage of adopting such non-dynamical 
simulations is that we can reach very large scales (we take 450h'^Mpc for the simulation box 
side), with a negligible computational cost, even when we average the results over several 
realizations. In a previous work (Amendola & Occhionero 1993), the bubble model has been 
compared to some observational constraints, such as the large-scale shape of 1 + ^(r) (here, 
^(r) is the galaxy spatial 2-point correlation function) and the upper limits on anisotropy 
and spectral distortion of the cosmic microwave background. In particular, it was found that 
power-law bubble models can fit the observed spatial correlation function if 

Rm ^ (p/130)-i-3 (2) 

(tested in the range p = 4 — 11) and if about 90% of the total volume is occupied by the 
bubbles. The latter requirement fixes the cut-off size Rm, i-e. the radius of the smallest 
bubbles. All the power-law models we test in this paper roughly obey the prescriptions above. 
If bubbles were generated before decoupling, Amendola & Occhionero (1993) also show that 
power law spectra pass the cosmic microwave background tests only for 6 < p < 11; roughly 
speaking, shallower spectra have too many large bubbles, while steeper spectra too many small 
bubbles. 

In Figure 1 we plot the resulting galaxy distribution for different bubble models, having 
both fixed and variable radii, according to the power-law spectrum of Eq.(|l]). Each panel 
reports the projection in the x — y plane of galaxies in a 15/i~^Mpc slice in the z direction. It 
is apparent the remarkable variety of structures generated by the bubble geometry, with rich 
clumps, filaments and underdense regions having sizes of some tens of Mpcs. Also the effect 
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of taking different model parameters is visible in this plot. Since each realization is generated 
with the same phase assignment, the persistence of several structures can be directly compared 
along the sequence of the model parameters. Fixed radius and steep spectra produces more 
isolated and concentrated galaxy aggregates, while shallower spectra turns into the presence 
of larger bubbles and, consequently, of more extended coherent structures. These differences 
will be quantified in more details in the remainder of this paper. 
The plan of the paper is as follows. 

After generating artificial Lick maps by projection of the spatial simulations, we analyze 
in Section 2 the resulting angular correlation functions and compare them with that observed 
for the APM sample (Maddox et al. 1990). Section 3 is devoted to the analysis of moments 
of cells counts (variance, skewness and kurtosis). The comparison with analogous results for 
observational data sets allows us to further restrict the parameter space for the permitted 
bubble models and to assess the reliability of this geometrical interpretation of large-scale 
clustering. In Section 4 we discuss our results and draw the main conclusions. 

2 Angular correlation 

As a first test of large-scale clustering, we investigate the projected (angular) correlation 
function implied by the bubble geometry. In order to implement this analysis, we generate 
artificial Lick galaxy catalogues, following the same procedure adopted by Coles ct al. (1993, 
hereafter CMPLMM) and Moscardini et al. (1993). The observer is placed at the center of 
the simulation box. Then, the original box is refiected, so to generate a three dimensional 
observational cone having width and depth adequate to reproduce the Lick map. Although 
the characteristic depth of the Lick map is ~ 210 h"^ Mpc, also very luminous galaxies at 
a much larger distance are included. Taking for the luminosity function the Schechter-like 
expression 



with the parameters given by Efstathiou, Ellis & Peterson (1988; a = —1.07 and M* = —19.68 
for the absolute magnitude associated to the characteristic luminosity L*), few galaxies up to 
~ 700 Mpc are included. Taking a simulation cube of 450 h"^ Mpc aside, we select galaxies 
up to the distance of 675/i~^Mpc. With this kind of construction the observational cone 
encompasses the range of galactic latitude > 45°. The procedure to assign luminosities to 
the galaxies is the same as described by CMPLMM. We suitably truncate ^{L) at the faint 
and bright tails, assign the apparent magnitude to a galaxy at a given distance by including 
X-correction and curvature effects. Therefore, we project only those galaxies brighter than 
the limiting apparent magnitude of the Lick map (m < 18.8). We tune the mean galaxy 
separation in the simulation box so to end up with a total number of projected galaxies of 
the same order of that included in the Lick map for > 45° (~ 316, 000). As observed by 
CMPLMM, although the observational cone includes galaxies belonging to different replica 
of the same simulation box, no substantial spatial periodicity should affect the projected 
clustering pattern, since each galaxy has a very small probability to be selected more than 
once. This problem is even less important in our case, being our simulation cube much larger 
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than that used by CMPLMM, which used three levels of replicated boxes. After generating 
the projected galaxy distribution, we collect them in counts with 10 x 10 arcmin cells, so to 
reproduce the observational setup of the Lick map. 

Simulating angular galaxy samples, like the Lick map, has been proved to be a reliable 
approach to compare data and simulations. In fact, the lack of redshift information in angular 
analysis is largely compensated by the extension of the available data sets. CMPLMM gen- 
erated artificial Lick maps from N-body simulations, based on assuming both Gaussian and 
skewed CDM initial conditions, in order to investigate the topology of the two-dimensional 
galaxy distribution. Moscardini et al. (1993) used the same approach to investigate the angu- 
lar correlation functions of such simulated samples. Borgani et al. (1993) extracted synthetic 
samples of galaxy clusters from the same Lick map simulations and applied several statistical 
tests to compare them with clusters selected from the real Lick map (Plionis, Barrow & Frenk 



Here, we analyze the angular correlation function, w{i}), of our Lick simulations and com- 
pare them with the APM data, as provided by Maddox et al. (1990). The adopted estimator 
of the angular two-point correlation function is 



Here, rii and rij are the galaxy counts in the i-th and j-th cell, respectively, placed at sepa- 
ration ■}}. The average is taken over all cell pairs with separation i?. Since our simulations of a 
bubbly Universe contain only geometry, we are not interested in the small-scale behaviour of 
w{'d), which is expected to be determined by non-linear gravitational dynamics. Considering 
only separation above 1?5, which corresponds to ~ 5 Mpc at the Lick depth (Maddox et 
al.l990), allows us to considerably reduce the computational time. We group the small cells 
to form counts in 1° x 1° cells, and w{{}) is evaluated according to eq.(|), being rii and rij the 
counts in such larger cells. 

The results of our analysis are shown in Figure 2, where we plot the angular correlation 
function for a list of bubble models, as compared to the APM data, suitably rescaled to the 
Lick depth. In Table I we report the relevant parameters for the models we are consider- 
ing. The upper left panel is for the fixed-radius model. It is apparent that the model with 
Rf = 12 Mpc displays a deficit of correlation strength at separations > 3°, while providing 
the correct small scale amplitude. Increasing the radius at Rf = lAh'^ Mpc increases w{'d), 
thus providing an excess of small-scale clustering. The upper right panel is for the variable 
radius models, with power index p = 8 and for two different values of the maximum bubble 
radius. It is apparent that allowing for a variation of the bubble radius sensibly improves 
the agreement. The shape of w{'d) is now quite well reproduced at all the relevant angular 
scales. In the lower left panel we check the effect of taking an even wider spectrum of bubble 
radii, by lowering the p value. Also in this case, the overall shape of w{-t}) is remarkably well 
reproduced. Finally, in the lower right panel we verify the effect of taking different amounts of 
unclustered galaxies within a fixed model. The plot shows that the w{{}) shape is left substan- 
tially unchanged apart from a variation of the small scale-amplitude. Therefore, increasing 
the percentage of unclustered galaxies amounts essentially to shift w{'d) downwards. This sug- 
gests that, once a model is shown to produce the correct correlation shape, its normalization 
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can be fixed a posteriori by choosing a suitable amount of random galaxies. It is however 
clear that this procedure is not completely arbitrary, since observational constraints on the 
percentage of weakly clustered galaxies in voids suggests that they should be ^ 10% of the 
whole population. Note that several models produces a w{i!}) slope at the smallest considered 
scales, which remarkably matches that for the APM sample. Since our simulations do not 
include gravitational dynamics, which is expected to be responsible for the small-scale w{'&) 
profile, this may suggest that the bubbly geometry can be rather adequate also to account for 
mildly non-linear clustering. 

The above results about angular correlation analysis indicate that the geometrical pattern 
provided by a random bubble distribution is able to provide a rather adequate amount of 
the large-scale power traced by the largest available projected galaxy samples. However, 
although several models are clearly ruled out, and others, such as that with fixed radii, are 
only marginally acceptable, the wli)) analysis seems not to be able to pick out one preferred 
model. For instance, looking at the plots of Figure 2, it is not clear whether a steep spectrum 
of bubble radii (with p = S) with rather small Rm is better than a having p ~ 5 and a larger 
Rm value. In order to further discriminate between different models, we apply in the following 
the higher-order statistics of count-in-cell moments to the three-dimensional distributions. 

3 Moments of cell-counts 

3.1 The method 

The analysis of moments of count in cells has been recently widely employed to analyze the 
observed galaxy distribution and it has been found to be a powerful statistical test of large- 
scale clustering (Efstathiou et al.l990; Saunders et al.l991; Bouchet, Davis & Strauss 1992; 
Coles & Frenk 1991; Gatzaiiaga 1992; Loveday et al.l992). The main advantage of this kind of 
analysis lies in the fact that it is relatively easy to implement and turns out to be particularly 
suited to investigate the statistics at large scales, where statistical noise severely limits the 
usage of correlation functions. 

In order to characterize the statistics of the moments of counts, let us introduce Sr^x) = 
Pi?(x)/p — 1 as the relative fiuctuations of the density field Pij(x), smoothed over a suitable 
angular scale Rhy a window function PFr(x). Therefore, 



being = J d'^xWji{'x.) the volume encompassed by the chosen window. At the scale R, the 
statistics of the fiuctuation field can be described by using the moments Hn{R) = {^r}- Let 
us consider 5r as a random variable with zero mean and let -P(5r) be its probability density 
function. Accordingly, its moment of order n is 




(5) 




(6) 



The corresponding generating function 




(7) 
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is defined in such a way that the /i„ moments represent the coefficients of its McLaurin 
expansion, 



oo 



(8) 

<f>=0 



An equivalent statistical description can also be given in terms of the cumulant generating 
function 



oo 

1^ # I I /I I I 

(9) 

<j>=0 



The cumulants k„ can be related to the moments /i„ by subsequently differentiating the 
corresponding generating functions. At the lowest orders, it is 

= /i2 ; = fJ-S ; K4 = - 3/^2 ■ (10) 

with more cumbersome relations holding at higher orders. In turn, the cumulant of order n 
represents the average value of the connected n-point correlation function C,n inside the window 
volume Vr: 

^-(^) = TTZvT / (nc?'a;,iy^(x,)Un(xi,...,x„). (11) 

Under the assumption that the 2-point function behaves like a pure power-law, ^(r) = (ro/r)"^, 
at the second order the variance is o"^ = ^2 = </2 where 

•^2 = 777^ / d""^! I f/'x2|xi-X2rWi(xi)iyi(x2) (12) 

[ViY Jvi Jvi 

arises after changing the variables in the integrals of eq.(|ll]), according to Xj :Ki/R. In this 
way, Vi represents the volume encompassed by the window of unit radius, Wii^), so that J2 
is a dimensionless coefficient, which only depends on the window proffie and on the slope 7 of 
the 2-point function. 

At higher orders, a popular model to describe the scaling of the moments is represented 
by the hierarchical expression 

Kn = Sn ^2 (13) 

where //„ = n — 1 and Sn are suitable coefficients, which depends on the moment order, other 
than on the shape of the window function. For instance, the hierarchical 3-point function 
reads ^3 = 3Q^|, so that eq. (pi]) gives % = 3(5(^3/^1)^2) where 



•^3 = TTTT^ / I nC?^a;il^i(Xj) I |xi-X2| ^|xi-X3| '^|X2-X3| 

(^1)^ Jvi L=i J 



(14) 



plays the same role as J2 for n = 3. On a theoretical ground the hierarchical recurrence relation 
is expected to be dynamically generated by strongly non-linear gravitational clustering ( e.g., 
Peebles 1980). In this regime, the closure of the BBGKY equations provides hierarchical 
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correlations with specified Sn coefficients (Fry 1984a; Hamilton 1988). Also in the mildly non- 
linear regime, second-order perturbative approaches to the evolution of gravitational clustering 
generates a hierarchical sequence of correlations ( e.g., Fry 1984b; Juskiewicz & Bouchet 1992), 
but with values of the coefficients Sn which are in general different from those arising from 
non- linear gravitational clustering (Lahav et al.l993). 

Suppose now to represent the continuous density field with a discrete point distribution, 
as it occurs for the galaxy distribution. In this case, precise relations exist to connect the 
moments of the continuous field, /i^, to those, m„, for the counts of the discrete realization. 
Under the minimal assumptions of Poissonian sampling of the continuous field, the discrete 
nature of the point distribution is accounted for by the usual change of variable (p —>■ e*^ — 1 in 
the functional dependence of the generators M(0) and K^cj)) ( e.g., Layzer 1956; White 1979). 
As a consequence, new terms appear when differentiating them, and eqs. (p!0|) become 

K2 = ^2 — 1/N (variance) ; 

^^3 = ""^3 ~ 3m2/iV — 1/iV^ (skewness) ; 

K4 = 7714 — 

67713/N - 3m^ + llma/iV - (kurtosis) (15) 

(Peebles 1980). Note, however, that the Poissonian representation to account for discreteness 
effects does not always provide a good prescription (Coles & Frenk 1991; Borgani et al.l993). 
This is the case, for instance, when the average count within is much less than unity, 
so that no informations about the underlying continuous statistics can be recovered. A fur- 
ther example is represented by the distribution of density peaks, which is far from being a 
Poissonian sampling of the background continuous field. 

In our analysis we take the window to be a cubic cell of size R. In this case, if we make 
a partition of the simulation box into A^";, cubic cells, the corresponding moments of counts 
reads 

rrin = — y^^-^ -, 16 

where A^ is the average cell count. Other authors considered different windows. Davis & 
Peebles (1983) analyzed the variance for the count of CfA galaxies inside "sharp" spheres, 
while Saunders et al.(1991) evaluated the variance and the skewness for the QDOT IRAS 
sample for counts within Gaussian spheres. In the following we will compare our analysis to 
analogous results from QDOT galaxies. Therefore, we need a suitable prescription to identify 
the size of the cubic cell to be associated to that of the Gaussian sphere. As for the variance, 
its value evaluated inside a Gaussian sphere of radius Rgs coincide with that evaluated in a 
cubic cell of side 

Rcc = {J2,cc/J2,gsf'^Rgs = 3.35 i?,, ; 7 = 1-8, (17) 

with obvious meaning of the indices, quite independently of 7. In general one expect a different 
relation for the scales connecting the skewness estimate within the two different windows, but 
we find numerically a relation very close to (|T^ within the MonteCarlo quadrature accuracy, 
again with a negligible dependence on 7, 

Rcc = {h,cc/Jz,gsf'^^Rgs = 3.25 i?^, ; 7 = 1-8, (18) 
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Note that, although the two prescriptions to identify the window sizes to be associated to a 
given moment value do depend on the moment order itself, in the case of interest we find only 
a rather weak dependence. 

We focus now on the results for a representative subset of the models, which are in fair 
agreement with the angular APM correlation function. The three fiducial models are F.14.3, 
8.38.1, and 5.66.1 (see Table I for the parameters relevant for these models). 

3.2 Results 

In Figure 3 we plot the scale-dependence of the variance for the three bubble models consid- 
ered, as compared with observational data. The error bars for our models are not reported here 
for clarity; they have been calculated as the ensemble scatter over ten different realizations, 
and are plotted on the following figures. Open circles are for our models, while filled circles 
are for the Stromlo-APM galaxies (Loveday et al.l992) and filled triangles for QDOT galaxies 
as analyzed by Efstathiou et al.(1990; hereinafter QDOT90). The agreement on large scales 
is fairly good, while a residual discrepancy is always detected at small scales. The reason for 
this could be that the variance for the APM galaxies has been evaluated in redshift space and 
then corrected to real space according to the relation a^edsh. sp. = fi^^^)'^reai sp.y where the 
correction factor 



holds at the scales of linear clustering and gives an equal amount of clustering amplification 
at all the scales (Kaiser 1987). In eq. ([l9|) , b is the linear biassing factor, which describe the 
difference between relative fluctuations of galaxies and underlying dark matter, 6gai = hdoM- 
Loveday et al.(1992) corrected their variance estimate by taking f{Q,b) = 1.4 at all the 
scales. However, this could underestimate the clustering at small scales, where, instead, virial 
motions are expected to suppress the redshift-space correlation ( e.g.. Fisher et al.l993). It 
is however not clear how large this effect should be in our case, since the bubble simulations 
does not provide galaxy peculiar velocities. In any case, we do not believe that this small- 
scale discrepancy represents a serious problem. A better agreement could be reached by 
allowing for a further 10% of unclustered component, which should decrease the values, 
without significantly changing the angular correlation (see Figure 2). However, here we are not 
interested in finding the model which best fits the observational data, rather we are analysing 
the ability of bubble geometries to reproduce the overall features of the galaxy distribution. 

In Figure 4 we report the relations of skewness fis and kurtosis /i4 versus the variance o"^. 
The models are compared with the corresponding relations found in the analysis of Gatzafiaga 
(1992) for CfA and SSRS galaxies (dashed lines) and of Bouchet, Davis & Strauss (1992) for 
the Strauss et al.(1990) IRAS 1.2Jy sample (dot-dashed lines, hereinafter IRAS 1.2). The 
scales analyzed in the two observed samples are limited to about 30 h~^Mpc for Gaztanaga 
(1992) and to about 70 /I'^Mpc for IRAS 1.2. We also plot the skewness -variance relation 
found by Saunders et al.(1991) in the QDOT survey with Gaussian windows (hereinafter 
QD0T91). The solid lines are the weighted best fit to our model's values, considering only 
the points for which o"^ < 1. The fitting parameters are reported in Table I; considering also 




(19) 
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the points with > 1 would result in values still inside the errors. We immediately note 
that, although the shallower power-law model 5.66.1 matches the APM angular correlation 



function, it is unable to follow the hierarchical relation (13) in the fully linear regime. 

Despite the apparently fair fit provided by the hierarchical relations. Figure 5 shows that 
we cannot exclude deviation from hierarchicity, due to the large error bars in both our model's 
and observational data. It appears indeed that the plots of skewness and kurtosis as functions 
of the variance may conceal even substantial departures from hierarchical scaling, depar- 
tures which, on the contrary, are best evidenced by plotting the coefficients S3 = ks/ct^ and 
5*4 = K^/a^ against the box scale. Again open circles are for our simulations while filled 
triangles are now for the QD0T91 analysis by Saunders et al.(1991), with errors propagated 
from those originally quoted. Dashed lines and dot-dashed lines denote the surprisingly small 
2a uncertainty level for the analysis of Gatzanaga (1992) and the la level for IRAS 1.2, re- 
spectively. We confirm with greater evidence that a shallow spectrum of bubble radii does not 
provide hierarchical correlations at large scales, and is totally at variance with the result at 
largest scale probed in QD0T91. On the contrary, steep power-law models and fixed radius 
models show a relative independence of S3 with the scale, at least within the statistical uncer- 
tainties, while the coefficient 5*4 remains completely undetermined at scales above 30/i~^Mpc. 
This does not allow a further discrimination among the models. 

Finally, in Figure 6 we display the results for ten simulations of the models 8.38.1 and 
F.14.3 in the skewness- variance plane. These are compared to the result from the QD0T91 
analysis of Saunders et al.(1991) with the largest Gaussian sphere, Rgs = 20/i"^Mpc. Adopting 
the relation Rcc = S.SRgg one gets Rcc = 66h~^Mpc, but since this scale identification depends 
on a model for ^2('") which, at these scales, is quite questionable, we plot the results for the 
ten simulations taking box sizes of 60 h~^Mpc (the rightward points) and 70 h~^Mpc (the 
leftward points). The satisfactory agreement once again stresses the ability of certain bubble 
models to account for the large scale clustering of galaxies. Notice, however, that the results 
in QD0T91 have not been corrected for the redshift distortion; applying a correction factor 
f{Q,b) > 1 would push down the variance by the same factor and the skewness by f{Q,by. 
Then the best agreement with our models would require somewhat less power on large scales, 
to be obtained, for instance, by increasing the unclustered component. 

4 Discussion and conclusions 

In this paper we realized a detailed clustering analysis of simulations of bubbly galaxy dis- 
tributions. Due to the purely geometrical character of our simulations, we do not include 
small-scale gravitational clustering. Rather, we concentrate on large (> 5/i~^Mpc) scales, 
where the clustering should be in the linear regime or weakly non-linear regime, so that it can 
be modelled by the geometry of the bubble distribution. In turn, the bubbly geometry can be 
thought to be either the result of a primordial phase transition or of quasi-linear gravitational 
dynamics, acting at intermediate and large scales. The advantage of this kind of approach is 
that we are able to generate very large simulation boxes of 450/i~^Mpc aside, which allow for 
a careful clustering analysis at large scales. 

As a first test, we generated artificial Lick maps, by projecting galaxies contained in two 
levels of replicated boxes (see CMPLMM and Moscardini et al.l993). The resulting Lick 
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simulations include galaxies up to a distance of 675/i~^Mpc and is limited to the range of 
galactic latitude > 45°. The possibility of reproducing the observational setup of this 
angular sample enables us to realize a comparison between the bubbly large-scale clustering 
and that of observational data sets. As a remarkable result, by comparing the resulting angular 
two-point correlation function, w{'d), with that of the APM sample (Maddox et al.l990), we 
find that its large-scale behaviour can be fairly reproduced by several models. Models with 
fixed radius, which are in any case less motivated, are found to have troubles. In fact, choosing 
Rm = 12/i~^Mpc, the resulting w{i)) amplitude agrees at t? ~ 2° with APM data, but fails at 
larger scales. Vice versa, for larger radii the agreement at 'd<^ 5° is sensibly improved, but it 
turns into an excess of clustering at smaller angles. Taking variable bubble radii introduces 
a modulation of the correlation amplitude, which succeeds at reproducing the observed w{i)) 
at all the relevant scales. However, we are not able to discriminate at this level between steep 
and shallow spectrum of radii. 

In order to further discriminate between different models, we also applied the statistics 
of moments of count in cells and compared the results with similar analyses on observational 
data sets (Saunders et al.l991, Bouchet, Davis & Strauss 1992, Gatzanaga 1992, Loveday 
et al.l992). As for the scale-dependence of the variance, we find that all the models, which 
satisfy the angular correlation test, are also able to account for the Stromlo-APM redshift 
data (Loveday et al.l992). This is not surprising, since both angular and spatial analyses 
of second order correlations are based on the same parent sample. More interesting results 
appear as we consider at higher correlation orders the scaling of skewness and kurtosis. In fact, 
the bubble geometry turns out to develop a rather well defined hierarchical scaling of higher 
order cumulants up to scales ^ 70/i~^Mpc, with the exception represented by the large-scale 
behaviour of the shallower spectra 4.78.1 and 5.66.1 (see Figures 4 and 5 and Table I). This 
result may appear rather surprising, since hierarchical scaling is expected to have a dynamical 
origin, both in the mildly and in the strongly non-linear regime of gravitational clustering, 
while no dynamics is present in our simulations. 

It is however possible to show that hierarchical scaling is a necessary consequence for any 
non-Gaussian distribution (like that implied by a bubbly geometry) in the weak coherence 
regime. To show this, let us assume to be at the scale r in a regime of weak clustering for 
the fluctuation field ^^(x). Let us also consider the coarse-grained field Sr{x.), smoothed at 
the scale R ^ r, that can be considered as the linear superposition of many independent 
processes ^^(x): 

Sni^) = - j2Sr{^^)■ (20) 

The corresponding moment of order n is 

l,^{R) = (5S(x)) = lEk^Ml. (21) 

Since 6r values are assumed to be nearly independent, the average of all the cross products 
in the multinomial expansion of eq. (|2T|) becomes negligible, so that fin{R) = P~*^"~^Vn(^)- 
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Therefore the same relation also holds between the cumulants and, since cr'^{R) = A*2(-R) = 
/X2(r)/p, it implies the hierarchical scaling 

Kr^iR) = 5„ [a(i?)]2("-^) . (22) 

with coefficients S'„ = Hn{r)/[o'{r)]'^^"~^^ fixed by the statistics at the reference scale r and 
independent of R. 

It is clear that this argument to explain the presence of hierarchical scaling crucially 
depends on the assumption of weak coherence. In most cases, this simply means that the 
clustering is well in the linear regime, (t|. ^ 1. This remarkably agrees with the results 
presented in Tabic I for models having different percentages of unclustered galaxies. In fact, 
as the unclustered component is increased, so to decrease the correlation amplitude, the 
reliability of the hierarchical scaling improves (see the values of 773 and 7)4 for the three models 
with p — 4.5). However, if the field (5(x) keeps coherent over large scales, e.g. due to the 
presence of intrinsic phase correlations, coherent structures can still be present even in the 
small variance regime and the hierarchical scaling does not show up. This is the case, for 
instance, when large-scale phase correlation is imprinted in the initial conditions or when 
dealing with primordial spectra providing a lot of large-scale power, which rapidly develop 
extended structures in the clustering pattern. This may explains, for instance, why in the 
analysis of N-body simulations by Bouchet & Hernquist (1992) it is found that, while the 
hierarchical relation is closely followed by the CDM model in the mildly non-linear regime, 
the evidence of hierarchical behaviour on large scales for the HDM model is much less clear. 
This is also the case for our 4.78.1 and 5.66.1 models, which do not develop hierarchical 
behaviours at large scales. In fact, taking a shallower spectrum of bubble radii allows the 
presence of very large bubbles, which generates a high degree of large scale coherence. 

In the introduction of this paper we addressed two questions concerning whether clustering 
models based on bubble geometry are adequate to account for the observed large-scale texture 
of the Universe and whether the resulting clustering measures are sensitive to the choice of 
the model parameters. Based on the results of our analysis we can give an affirmative answer 
to both questions. It is however clear that our findings demand for a dynamical interpretation 
about the origin of a bubbly galaxy distribution. In particular one may ask whether they are 
imprinted as initial conditions in the primordial fiuctuations or develop through gravitational 
evolution. In this respect, one could be tempted to suggest that, while bubble generation 
through primordial phase transitions should give rise to a narrow spectrum of radii (which is 
also more likely to pass the cosmic microwave background constraints), gravitational evolution 
of a rather flat fluctuation spectrum should produce nearly spherical voids of sizes comprised 
in a large range. It is however clear that, before rush to definite conclusions on this point 
one should wait, from the observational side, for more extended data sets, which were able to 
encompass the largest scale involved by coherent structures, and, from the theoretical side, for 
more rigorous dynamical descriptions like those provided by running large N-body simulations. 

Acknowledgements. We wish to thank Dr. Steve Maddox for providing us with the digitized 
version of the APM correlation data. 
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Figure captions 

Figure 1. Slices of the bubbly galaxy distribution for different models (see Table I for the 
parameters of each model). Each slice is 15h~^Mpc thick and all the realizations correspond 
to the same choice of initial phase-assignment. 

Figure 2. The angular correlation function, w{i)), for simulations of Lick map, obtained by 
projecting the three-dimensional bubble simulations. Open circles refer to the APM corre- 
lation, rescaled to the depth of the Lick map, as provided by Maddox et al.(1990). Filled 
symbols are for our simulations. They are plotted only for t? > 1.5°, which correspond to the 
physical scales of linear clustering at the depth of the Lick map. 

Figure 3. The scale dependence of the variance for counts in cubic cells. Open circles are 
for the bubble simulations, which are obtained after averaging over ten realizations. Filled 
triangles are for the IRAS QDOT analysis by Efstathiou et al.(1990) and filled circles are for 
Stromlo-APM galaxies (Loveday et al.l992). For reason of clarity, we do not plot errorbars 
for our data. 

Figure 4. The variance-skewness (left panels) and the kurtosis-skewness (right panels) rela- 
tions. Open circles refer to our simulations, with errorbars denoting the r.m.s. scatter over an 
ensemble of ten realizations. Heavy bars for the skewness are the QDOT data by Saunders 
et al.(1991). The dot-dashed lines are the best fits to the analysis of the IRAS 1.2Jy survey 
by Bouchet, Davis & Strauss (1992) and dashed lines refer to the analysis of CfA and SSRS 
samples by Gatzafiaga (1992). The solid lines are the best fits to our data, according to the 
parameters reported in Table I. 

Figure 5. The scale dependence of the (left panels) and 6*4 (right panels) hierarchical 
coefficients of eg. ([T3|) . Open circles are for our models, while dot-dashed and dashed lines 
delimit the uncertainty bands in the IRAS 1.2Jy analysis by Bouchet, Davis & Strauss (1992) 
and Gatzafiaga (1992), respectively. Filled triangles for 5*3 refer to the Saunders et al.(1991) 
analysis of QDOT data. 

Figure 6. The variance-skewness relation for the QDOT data at the largest considered scale, 
corresponding to the size of 66/;,~^Mpc for a cubic box. Each open circle refers to the result 
for one realization of the corresponding bubble model. Leftmost circles are for a box-size of 
70/i~^Mpc, while rightmost ones are for the 60/i~^Mpc scale. 
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